Testing the OPE Wilson coefficient for A 2 from lattice QCD with a dynamical charm 
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Gluon and ghost propagators data, obtained in Landau gauge from lattice simulations with two 
light and two heavy dynamical quark flavours (AT/=2+l+l), are described here with a running 
formula including a four-loop perturbative expression and a nonperturbative OPE correction dom- 
inated by the local operator A 2 . The Wilson coefficients and their variation as a function of the 
coupling constant are extracted from the numerical data and compared with the theoretical ex- 
pressions that, after being properly renormalized, are known at 0(a 4 ). As also A^g is rather well 
known for Nf =2+1+1, this allows for a precise consistency test of the OPE approach in the joint 
description of different observables. 

PACS numbers: 12.38. Aw, 12.38.Lg 

I. INTRODUCTION 

The Green functions of a general quantum field theory contain in principle all the information of the theory. 
Consequently, since a few decades, considerable research effort is being devoted to construct as reliable as possible 
estimates of these Green functions in QCD. The main motivation for this has been to get some understanding 
about the infrared properties of QCD (see pQ for a recent review), but they have been also aimed at extracting 
low-energy parameters as Aqcd from lattice QCD computations (see, for instance, [2-8 ). In particular, gluon and 
ghost propagators have been recently computed, within the ETM collaboration, from lattice gauge fields including 
Nf— 2+1+1, two light degenerate and two heavy, dynamical quark flavours (see [HI HO] for details of the simulations 
set-up) with a twisted- mass fermion action [TTJ [T^] . Those propagators were first exploited for the calculation of the 
MS Strong running coupling via the Taylor-scheme coupling in Landau gauge [Ml [14] and, more recently, with pure 
Yang-Mills (no flavours) and two-light-flavours propagators, all together studied in order to investigate the effects on 
them from dynamical quark flavours |15) . The running coupling estimated through these propagators with Nf =2+1+1 
at r and Zq mass scales has been proven to agree pretty well with experimental results from r decays |16j and " world 
average" from PDG [17]. In the analysis of refs. [HI [14], the invocation of nonperturbative OPE corrections is 
unequivocally needed to account properly for the lattice data of the coupling in the Taylor scheme. 

The OPE approach allows for the expansion of the matrix element of any non-local operator in terms of local 
operators, that can be conveniently organized in a hierarchy by their momentum dimensions. Then, the terms 
of that OPE expansion for any Green function provide with nonperturbative corrections that, after the sum rules 
factorization jTB] US], appear coded as a coefficient to be computed in perturbation (Wilson coefficient) and the 
nonperturbative condensate of a local operator. Furthermore, a consistent OPE picture implies that, with the proper 
renormalization of the local operators, their condensates should be "universal" and take the same value in the OPE 
expansion for any different Green function. In refs. [HI [14] (also in refs. [20] for Nf=0 and Nf —2 quark flavours), 
the Taylor-scheme running coupling is confronted to its OPE prediction to leave us with a best fit of both Aqcd and 
the dominant nonperturbative condensate, i.e. that of the local operator A 2 which has been elsewhere very much 
studied {e.g., see [21TI27] ). It is worth to emphasize that, in constructing the Taylor-scheme coupling, the bare lattice 
propagators appear combined such that the cut-off dependence must disappear and any correction depending on the 
lattice spacing hence will vanish when approaching the continuum limit 6]. This makes easier the confrontation of 
lattice data and continuum predictions but prevents the check for the universality of condensates, as gluon and ghost 
propagators are not to be independently analyzed. On the other hand, concerning the universality of condensates, 
previous studies exploiting quenched (Nf=0) [7ll51 l2"Tll2"8] - l50"] lattice gauge configurations have shown a fair consistency 
for the values of this dimension-two gluon condensate extracted from the dominant nonperturbative corrections for 
gluon, ghost and quark propagators and the three-gluon vertex, in Landau gauge. The same have been found when 
comparing the condensate extracted from lattice quark propagator and the Taylor coupling with Nf=2 dynamical 
flavours [201 EI] ■ Still, despite that consistency and the pretty convincing estimate of Aq CD obtained with the good 
agreement of the OPE formula and the lattice results for the Taylor coupling with Nf =2+1+1, the nature itself of 
the unavoidable nonperturbative corrections and the involved condensates is a controversial matter (there was, for 
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instance, some recent controversial work about the notion of condensates that can be followed in refs. [3"2"rf51] and 
references therein). 

In the present note, we will make one more step with respect to refs. |13l 114) and will face any controversy about 
the nonperturbative OPE approach only by following that " the proof of the pudding is in the eating" . Then, our goal 
here will be to test again the validity of the OPE approach, now describing gluon and ghost propagators computed in 
Landau gauge from realistic lattice simulations including non-dcgcncratcd strange and charm quarks (Nf=2+1+1). 
The running with momenta for both propagators is well known at the four-loop level from perturbation theory. As 
one deals with a very realistic simulation of QCD at the energy scale for the r physics, Aqcd will be taken as an input 
from experiments with r decays. Therefore, the deviation with respect to the perturbative running for the lattice 
propagators will be accurately isolated and identified to the leading OPE power correction. This will allow for a very 
precise comparison that will be sensitive to the running with momenta of the Wilson coefficient for the local operator 
A 2 , which is known at 0{a A ) order. Thus, the deviation from the perturbation theory for gluon and ghost lattice 
propagators will be accurately accommodated within a nonperturbative OPE correction, including the running of the 
Wilson coefficient, by fitting a value of the condensate for each propagator and verifying that they both agree. This 
is a double and very demanding check. 

The paper is organized as follows: the nonperturbative running formulae for the propagators, after properly defining 
the renormalization prescription, are derived in section [n] the results are compared with their lattice estimates in 
section |TTT| the Taylor coupling is again obtained from bare gluon and ghost propagators and, for the sake of consistency, 
compared with the nonperturbative prediction in section |TV] and the conclusions are finally presented in section [V] 



II. THE WILSON COEFFICIENTS AT O (a 4 ) 

The purpose of this section is to derive the running of the propagators including an OPE contribution dominated 
by the local dimension- two operator A 2 . In particular, we focus on the precise definition of the renormalization 
prescription and in obtaining the Wilson coefficients at C(a 4 ) order. We will first compute the OPE-based equations 
for both gluon and ghost propagators and then combine them to derive the expression describing the nonperturbative 
running for the Taylor coupling. 



A. Two-point functions 

For the sake of generality, let us invoke a bare scalar form factor, T{p , A), defined from any gauge-dependent 
two-point Green function, where A is a regularization parameter and p 2 the momentum scale. We can formally write 
the following OPE expansion: 

t(p 2 A) = r pcrt (V a) (i + rA ( p2,A ^ ^ A + ... | m 

L[p,A) y P ,YV| ^ + rport (p2 ) A)4(JV|-l)p2 J 

which is taken to be dominated in Landau gauge by the local operator A 2 and where the purely perturbative con- 
tribution 1 is explicitly factorized. The bracket in r.h.s. of Eq. (nl) must be finite for a renormalizable theory, as the 
same multiplicative renormalization constant is required to remove on a same footing all the the terms in the OPE 
expansion. Furthermore, one can apply the dimensional regularization procedure (this implies A = 2/(4 — d) — £ _1 ) 
and then, for dimensional reasons, the term proportional to 1/p 2 cannot produce an additional pole at d = 4 in a 
renormalizable theory. Therefore, the poles for the obviously divergent quantity inside the bracket 



ppert (p2 j£ -l) 

and those coming from the condensate (A 2 ) E -i kill each other. In other words, when dealing with the UV divergencies, 
one can renormalize the ratio of Eq. ^ as follows 



1 Using the path integral jargon, it corresponds with the contribution to the Green function from field configurations with a null vacuum 
expectation value of the local operators in the OPE expansion. 
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where the dimension-two local operator is renormalized by Zjp such that 

(A 2 ) R ^ = Z-J{n\e- 1 ) (A 2 ) e -i (4) 

and we prescribe the standard MS renormalization recipe that only 1/e-poles and accompanying factors log(47r) are 
to be dropped away. Eq. |s|'s r.h.s. stands for the definition of c^ VMS , in the "Wilson Ms" scheme 2 (Ms) [3Tj . after 
the introduction of 

MOM / 2 2n _ rPClt (P 2 ' £ *) ( cr\ 



which allows to recast Eq. ([TJ) that now will read 



„WMS 



Thus, one can immediately notice that the renormalization of the Wilson coefficient inside the bracket and that of 
the two-point Green function outside are independent processes, irrespective with each other. In particular, one can 
also prescribe the standard MS renormalization scheme for the two-point Green function and will obtain the following 
helpful relation: 



(g,a(M 2 )) 
(^«(M 2 )) 



2 s / 2 . C WMS 

,a(M 2 )W S feXM 2 ))^— 4== (rl (7) 



fi 2 ' J \n 2 " 7 c mo \ j (p. 

that will allow for a direct link of c^ MS to the standard MS Wilson coefficient computed at O (a 4 ) -order in [35]. For 
the purposes we pursue, the running with momentum, p, for the OPE nonperturbative correction inside the bracket 
of Eq. (D can be derived by solving the following renormalization group (RG) equation 

"7^ S (a(M 2 )) + + P («(M 2 )) #4 Y { = , (8) 

' aiog/i 2 yy ^ " da J r¥ OM (K,a(n 2 )) 



1 o 



with 

the anomalous dimension for the dimension-two local operator A 2 that can be found at the four-loop order in ref. [36j . 
This RG equation is obtained just by applying the logarithm derivative on the renormalization momentum, /i, to the 
two hand sides of Eq. The boundary condition for the solution, 

cf^^ajp 2 )) c^(l,a(p 2 )) 

<g™(l,a(?)) c^(l,a(p 2 )) ' 1 ' 

is provided by Eq. |t| and can be borrowed from refs. [35 , where c 2 vls (l, a(q 2 )) is given by taking q 2 — in 

Eq.(7)/Eq.(9) for the gluon/ghost case, and [37], where c^ IS (l,a) for both ghost and gluon two-point Green functions 
is given in appendix C. 

Finally, we only need to invoke multiplicative renormalizability to extend the validity of Eq. (|6| to two-point Green 
functions obtained by applying any other regularization scheme. In particular, for lattice regularization, we will have 

fr? \ ( c™f4,«(p 2 ) N ) (A 2 ) 




2 This scheme is defined for the Wilson coefficient in the OPE expansion by imposing that the local operator is renormalized in MS, while 
the expanded Green function is in MOM scheme. 
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where a is the lattice spacing and z(/i 2 , a *) is an overall factor that, according to Eq. coincide with T peTt (fi 2 , a 1 ). 
In order to solve Eq. (JsJ , one can write 

C(a(|U ),a(p )) (12) 



„MOM / P 2 
c 



(g,a( M 2 )) c^(l,a(p 2 )) 



and solve for £(a(/z 2 ), a{p 2 )) the same Eq. ([8]) with the boundary condition ((a(p 2 ), a(p 2 )) — 1. Thus, we will obtain 
at the order 0(a 4 ) 



„WMS / £ 



(S.«(M 2 )) 



f OM (j,a(M 2 )) V«(M 2 ) 

(1 + 1.31808 a{p 2 ) + 1.51845 a 2 (p 2 ) + 3.04851 a 3 (p 2 )) T = G 



qY \ loo 

' 1 127ra(/i 2 ) (1 - 0.54993 a(/i 2 ) - 0.14352 a 2 (/x 2 ) - 0.07339 a 3 (/x 2 )) 



1 + 0.89698 a{p 2 ) + 0.84882 a 2 (p 2 ) + 0.87056 a 3 (p 2 )) T = F ' ^ 



where G and _F stand, respectively, for the gluon and ghost dressing functions and we take a = olt, the Taylor 
coupling, for the expansion. 



B. The Taylor coupling 



As the Taylor coupling is directly obtained from the lattice in terms of gluon and ghost dressing functions, up to 
lattice artefacts corrections, it follows jS], 



a T {p ) 



1 - ^IfV-"- 1 )^^ 1 ) 

47T 



Then, one can combine eqs. ( 11|13 ) and be left with 



(14) 



a T (p 2 



„pcrt / 2 



(p 2 ) i + 9i?K crt (f 2 ) J < rt (M 2 )) 



a(p 2 )\^ g 2 {„ 2 ){A 2 )^ 



a(/i 2 



4 (TV 2 - l)p 2 



(15) 



where 



R (a, a ) = (1 + 1.03735 a + 1.07203 a 2 + 1.59654 a 3 ) (l - 0.54993 a - 0.14352 a 2 G 



0.07339 al) 



(16) 



following result, 



In obtaining Eq. (15 1, the purely perturbative contribution, a^. clt , appears factored out of the bracket because of the 



F, 



MOM / P_ 



,a(n 2 ) 



r<MOM 
^0 



,a(/^ 2 ) 



nort 



(P 2 ) 



pert 



(m 2 ) ' 



(17) 



where c^ OM , defined in Eq. is here dubbed G^ OM for the gluon case (r = G) and F MOM for the ghost (T = F). 
Eq. (17) can be straightforwardly derived from the nonperturbative definition of the Taylor coupling by Eq. (fl4| but 



also proved at any order in perturbation. For the latter, one only needs the gluon and ghost anomalous dimension 
in MOM scheme and beta function in Taylor scheme to compute in perturbation F^ OM , Q^ LOM and a^ ert and then 
verify Eq. (|17|) order by order. 



The reader might have noted that coefficients in Eq. ( 16 ) appear to be slightly different with respect to those in 

in this 



Eq. (6) of rcf. [13 . This discrepancy is explained by the different boundary conditions used to solve Eq 
paper, Eq. (10), and in ref. [13], where the denominator of Eq. (10)'s r.h.s. is approximated by 1. As Eq. (6) of ref. [T3J 
is applied to estimate the strong coupling in ref. [T3], replacing Eq. (6) of [T3J by Eq. (16) would affect the strong 



coupling estimates, although with an almost negligible impact: the systematic deviation for A^i- 4 would increase it 



roughly by 1 %, i.e. 328(18) MeV instead of 324(17) MeV. 



MS 
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III. GLUON AND GHOST PROPAGATORS 



The analysis to be done now consists in applying Eq. (Ill with ( 12|13 ) to describe the lattice gluon and ghost 
propagators that have been computed for the running coupling study of ref. [14] and used for the study of dynamical 
quark flavour effects in ref. 15J, 



F ab (q 2 ) = \, (^expfa. (*-»)] (M" 1 ) 



q»qv\ G(q 2 ) 



(18) 



where A® is the gauge field and M = d^D 1 ^ is the Fadeev-Popov operator. The reader interested in details about 
the computation are referred to section II and III of ref. |15j . while here we will just analyse the results for gluon and 
ghost dressing functions, defined by Eq. (18) and obtained from the lattice simulations in p3] and described in Tab.JTj 
after the 0(4)-breaking lattice artefacts have been cured by the so-called i?(4)-extrapolation procedure [H[35]. These 
estimates for bare gluon and ghost dressing functions appear plotted in Fig. [I] 



p 




am 




a/is 


{L/af x T/a 


confs. 


1.90 


0.1632700 


0.0040 


0.150 


0.1900 


32 3 x 64 


50 


1.95 


0.1612400 


0.0035 


0.135 


0.1700 


32 3 x 64 


50 


2.10 


0.1563570 


0.0020 


0.120 


0.1385 


48 3 x 96 


100 



TABLE I: Lattice set-up parameters for the ensembles we used in this paper: K cr it is the critical value for the standard hopping 
parameter for the bare untwisted mass; /i; stands for the twisted mass for the two degenerated light quarks, while \i a and fis 
define the heavy quarks twisted masses; the last column indicates the number of gauge field configurations exploited. This 
implies that the strange quark mass is roughly set to 95 MeV and the charm one to 1.51 GeV (in MS at 2 GeV), while degenerate 
light quark masses range from 20 to 50 MeV (The lightest pseudoscalar masses approximately range from 270 to 510 MeV). 



A. Removing 0(4)-invariant lattice artefacts 



The gluon and ghost propagator lattice data shown in the upper plots of Fig. [Tjare still affected by 0(4)-invariant 
lattice artefacts. To treat these artefacts, we will apply the strategy that, in the following, will be described for the 
gluon propagator case. 



As discussed in ref. [T31 HH EI] for the cases of Taylor coupling and the vector quark propagator renormalization 
constant, one can properly parametrize the remaining 0(4)-invariant lattice artefacts as follows: 



Z 3 Lat V,0 = G{p\a- 1 ) + c a2p2 a 2 p 2 



(19) 



where Z\ att is the gluon propagator renormalization constant in MOM scheme which corresponds to the lattice gluon 
dressing function. In r.h.s., the artefacts- free gluon dressing is noted as G and, in the appropriate range of momenta, 
is supposed to be given by Eqs.([lljl3|). Then, Eq. Q can be recast to read 



Z Latt(p2 )a -l) 
GNP(V, a ( M 2) 



z(/i 2 ,a l ) + c a 



2p2 



2 2 

a p 



GNP(^,a( M 2 ) 



(20) 



with 



„WMS 



(A 2 ) 



„MOM / P 
L 



2] \ 4(iV c 2 -l)p 2 



+ 



(21) 



where 



WMS / JVIOM 



is given by Eq. (13|'s r.h.s. for T = G and where Gg , defined by Eq. (|5|'s r.h.s. again 



with r = G, can be computed at the four-loop order with the help of the perturbative expansions for the gluon 
propagator MOM anomalous dimension and Taylor scheme beta function [37] . Provided that G NP describes properly 



() 




FIG. 1: Bare gluon and ghost propagators as they result from lattice computations before (top) and after being cured of 
0(4)-invariant lattice artefacts and properly rescaled to appear superimposed (bottom). The blue solid line correspond to the 
best fit with the OPE formula. 



the nonperturbative running of the gluon dressing over a certain range of momenta, Eq. ( 20 ) can be used to make 
a linear fit of data within this range 3 for the three exploited data sets simulated at three different bare couplings, 
5o( a_1 ) = This will allow for the determination of the overall factor, z(n, a -1 ), for each data set and the 

coefficient c a 2 P 2, thought to be roughly the same 4 for any /? (see Fig. [2] and tab. [n]). Once c a 2 P 2 known, one can 
remove the 0(4)-invariant artefacts away from the lattice data for all available momenta, by applying Eq. (|19[ ), and 
be left with the continuum nonperturbative running for the gluon dressing; while the overall factors z(fj,, a _i ) can be 
used to rescale the data obtained at different j3 such that they all will follow the same curve (see down plots of fig. [I]) . 

As G NP depends on A^jg and the gluon condensate, g 2 (A 2 ), their values have to be known prior to the removal of 
0(4)-invariant artefacts. In practice, we will take Aj^g- to be known from r decays (a^g(m^) = 0.334(14) [IB], i.e. 
Ajvis = 311 MeV for Nf = 4) and will search for a value of g 2 (A 2 ) such that one gets the best matching of the rescaled 
gluon dressings from the three data sets. The criterion to define the best matching, that is explained in appendix, 
implies the minimization of certain parameter, A p , computed with the square of the area between polynomial fits of 



3 At lower momenta, the range is bounded by p = 4.5 GeV, below which the lattice data deviate from the nonperturbative prediction only 
including the OPE leading power correction. At large momenta, the upper bound is provided by a 2 p 2 ~ 3.5, above where the impact of 
higher-order lattice artefacts cannot be neglected any longer and results, for instance, in deviations from the linear behavior suggested 
by Eq. |20) . 

4 By dimensional arguments, the coefficient c a 2 p 2 might only depend logarithmically on the lattice spacing a(/9). Similar corrections, 
there inspired on the lattice perturbative expansions for the propagators, are applied in ref. |39) . 
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gluon 


ghost 


g 2 {A 2 ) [GeV 2 ] 


4.7 


3.2 


P 


2.10 


1.95 


1.90 


2.10 


1.95 


1.90 


2(^ 2 ,a _1 ) 


0.6495(5) 


0.6068(26) 


0.5971(30) 


1.215(6) 


1.229(20) 


1.193(31) 


Ca2p2 


-0.0097(2) 


-0.0071(8) 


-0.0080(8) 


0.0002(19) 


-0.0012(44) 


0.0011(69) 


g 2 (A 2 ) [GeV 2 ] 


4.1 


4.1 




0.6531(5) 


0.6106(26) 


0.6014(30) 


1.206(6) 


1.218(19) 


1.182(30) 


C a 2p2 


-0.0104(2) 


-0.0077(8) 


-0.0086(8) 


0.0022(19) 


0.0004(44) 


0.003(7) 



TABLE II: The overall factors, z(^i 2 ,a _1 ), and the coefficients c a 2 P 2, obtained as explained in the subsection III A for optimal 
values of g 2 (A 2 }. As discussed there, c a 2 P 2 is thought not to depend very much on f3 and so is fairly obtained from the results 
shown here. The renormalization point, fi, for the overall factor is fixed at 10 GeV. Statistical errors have been comput ed by 
applying the Jackknife method. The errors for the determination of the optimal values of g 2 {A 2 ) are also discussed in 
and can be found in next Tab. Mil 



III A 




FIG. 2: (Left) The parameter A p , related to the area between polynomial fits of the dressings for each two data sets (see 
appendix), in terms of the value for the gluon condensate, g 2 {A 2 ). The black solid line corresponds to the gluon, the black 
dotted to the ghost case and the red solid one to the sum of both. The areas shadowed in blue and brow n co rrespond to the 



1-ct allowed regions according to refs. Q3] and P3], respectively. (Right) Linear fits, as suggested by Eq. ( 20 1 , of Z^ tl /G NP 
(gluon) and Z^ a ' tt /G NF (ghost) in terms of a 2 p 2 /G NP , where g 2 (A 2 } takes the total minimum shown in the left plot, for the 
three lattice data sets. 

Thus, after applying the above described procedure to sweep the 0(4)-invariant artefacts away from data, we obtain 
the minima for A p in terms of g 2 (A 2 ) and the results for c a 2 P 2 and z(/i 2 , a -1 ) shown in Tab, pi} In the left plot of Fig. [2j 
the minima of A p for both gluon and ghost dressing functions are shown to take place at values for g 2 (A 2 ) that are 
fairly compatible with the estimate of [2]: 3.8(1.0); and with that of [13J : 4.5(4). For the sake of a proper comparison, 
one needs to estimate the uncertainty in determining the best value of g 2 (A 2 ) by the minimization of A p . To this 
purpose, the Jackknife procedure can be applied and gives for the gluon case the following result: g 2 (A 2 ) = 4.7(1.6) 
GeV 2 . In the ghost case, as the coefficient c a i P 2 is found to be negligible (see Tab. [n]), one can futhermore take 
c a 2p2 = and so avoid very large errors due to strong correlations between this coefficient a nd b oth z(/^ 2 ,a _1 ) and 
g 2 (A 2 }. Thus, we will obtain g 2 (A 2 ) = 3.1(1.1) GeV 2 . Those results appear collected in Tab."~ 
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Finally, one can impose the condensates to be the same for both gluon and ghost propagators and a total minimum 
for the sum of gluon and ghost A p 's will be found for g 2 (A 2 } =4.1 GeV 2 , that appears to be again in pretty good 
agreement with the previous results in refs. (T3l [14]. In the following, we will use the artefacts-free gluon and ghost 
dressing functions, obtained with this last value for the gluon condensate, shown in Fig. [I] 





gluon 


ghost 


a T H3] 


a T [T31 


g 2 (A 2 ) [GeV 2 ] 


4.7(1.6) 


3.1(1.1) 


4.5(4) 


3.8(1.0) 



TABLE III: Results for g 2 {A 2 } independently obtained from gluon and ghost propagators in this paper and from the Taylor 
coupling analysis in refs. [131 114] . 



B. Testing the OPE Wilson coefficients from data 



At this point we are in a position to extract directly from the numerical data the running induced by Wilson 
coefficients, obtaining thus their behaviour as a function of the coupling constant and making comparisons with the 
theoretical expectations, derived in sec. [TTJ The ability of the OPE to describe consistently and accurately different 
Green functions can then be tested. Such is a smoking gun for the reliability of the OPE nonperturbative approach. 
To this goal, Eq. (21 1 can be recast as follows, 



/-.MOM I P I 2 



r MOM I P 
U 



,a(/z 2 ) 



„WMS 



(A 2 



„MOM 

-0 



(g.«(M 2 )) 



A{N 2 



(22) 



where the l.h.s. is to be computed from the artefacts-free dressing function, G, defined by Eq. (19 1 and obtained, 



as explained in the previous subsection, by requiring the optimal matching of data from the three lattice data sets 
at different /?. On Eq. (22|'s r.h.s., the running wi th momenta can be easily expressed as a perturbative series in 
terms of the Taylor coupling (see, for instance, Eq. (13l), being proportional to the gluon condensate, g 2 (A 2 ), which 
has been also obtained in the previous section through the best-matching criterion. Then, one can plot Eq. |22l)'s 



l.h.s. in terms of ar and compare with r.h.s. 's prediction. This can be seen in Fig. [3] for the gluon (top) and ghost 
(bottom) cases. In view of this, we can conclude that: (i) The nonperturbative corrections are totally unavoidable, 
as data for Eq. (22)'s l.h.s. are clearly not compatible with zero; (ii) the Wilson coefficient derived at the order 
0(a 4 ) in sec. Ogives a high-quality fit of data (blue solid lines in both plots); (iii) had we not included the Wilson 



coefficient running (i.e., dropped c™ M 7 c away from Eq. (|22[)'s r.h.s.), or just included its leading logarithm 



approximation, a flatter slope would be obtained not describing properly the data (blue dotted and green solid lines). 
It is worthwhile to emphasize that when applying a simple linear fit of data within the appropriate range of momenta 
(i.e., 0.250 < o?t < 0.295), one obtains the orange dotted line shown in the right plots, the slope of which agrees 
pretty well with the local one resulting from the full r.h.s. of Eq. (22 1, in both gluon and ghost cases. 



IV. RECONSTRUCTING THE TAYLOR COUPLING 



After applying the procedure to sweep the (9(4)-invariant lattice artefacts away, we are left with the artefacts-free 
gluon and ghost dressing functions. One can then apply the definition, Eq. (14), to get an estimate of the Taylor 
coupling that is supposed not to be affected by artefacts at the order 0(a 2 p 2 ). In ref. [13], we proceeded the other 
way around: we first computed the Taylor coupling from lattice data through Eq. ( 14 ) and only then applied both 
the -ff (4) extrapolation and 0(4)-invariants artefacts removal to cure hypercubic artefacts at the desired order. Thus, 
re-obtaining the Taylor coupling, now after dropping the lattice artefacts away from the dressings, must be a strong 
consistency check for the artefacts treatment and for the whole analysis. As the artefacts-free dressing functions have 
been rescaled such that data for f3 — 1.90 and (3 = 1.95 result superimposed onto those for /3 = 2.10 (see Fig. [lj, we 
plug these rescaled dressings into Eq. ( 14 1 and take the bare coupling, g 2 , for (3 — 2.1. The results are shown in Fig. [4j 
where a fit with Eq. (15), both and g 2 (A 2 ) taken as free parameters, is also plotted. One obtains for the fitted 
values (including jackknife errors): A^ = 320(10) MeV and g 2 (A 2 ) = 3.9(3) GeV 2 ; strikingly in agreement with the 



results from the direct analysis of gluon and ghost dressings, in previous section, and with those from refs. [13\ I14|. 

In the following, we will emphasize that obtaining compatible values of A^jg and g 2 (A 2 ) from both a fit of Eq. ( 15 ) 
to the Taylor-coupling lattice data and fits of Eq. ^ to gluon and ghost dressings is neither trivial nor tautological. 
In the case of the renormalized Taylor coupling, once the artefacts been properly removed, not only the running with 
momenta is determined by the fitted parameters but also the size itself of the coupling at any momentum. To make 
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FIG. 3: Scrutinizing the Green functions data to isolate the running for the Wilson coefficients. The left plots shows the 
difference of the rescaled lattice data for the gluon (top) and ghost (bottom) dressing functions and their perturbative estimate, 
multiplied by the square of momentum, plotted versus the perturbative Taylor coupling at the same momentum. For pertur- 
bative evaluations, A^jg is taken to be known from r decays. Blue dotted and green solid lines respectively correspond to the 
leading logarithm approximation for Wilson coefficient and to its suppression. The right plots only show the region, from 4.5 
to 6.3 GeV where the behavior predicted by the OPE analysis in sec. [IT] is observed. The dotted orange line shows a linear fit 
of data that results to agree very well with the predicted running. 



it clear, Eq. ([6]) for the gluon (r = G) and ghost (r = F) dressings can be aso plugged into Eq. ( 14 ) and one gets 



/ 2\ 9o( a ) 2 11 -1\ / 2 — 1\ ( ciMOM ( P 

"r(P ) = 4?r z gh , a )z g i{fi,a ) \F I- 



«(M 2 ) 



r,a0 2 ) 



6 4«0* 2 ' a l ) 



z g i (/x 2 ,a x ) 



Att/3 



pert 



(m 2 



a 



pert 



(p 2 ) i + 



(23) 



where z g i (z g h) is the overall factor for the gluon (ghost) dressing function determined by fitting the 0(4)-invariant 
artefacts with Eq. (201 and reported in Tab. [TTj The brackets with dots inside trivially correspond to the same bracket 
accounting for the gluon condensate OPE correction in Eq. (151 and, in obtaining the second-line r.h.s., the latice 

As Eq. 



bare coupling definition, 



6//3, and Eq. (17) are also applied. As Eq. (Ill]) provides with a precise description 
from r decays and the gluon condensate reported in Tab. ITJ over 



for the gluon and ghost lattice dressings with A^jg 

4.5 < p < 6.3 GeV, Eq. (23l's r.h.s. necessarily provides with the momentum behaviour within the same window and 
for the same parameters. However, Eq. (fl5|) will be only re-obtained if 



pert 



(M 2 



6 



4h (/A a X ) z gl OA" *) 



(24) 
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FIG. 4: The Taylor coupling computed with the artefact-free lattice dressing functions and with Eq. ( 15 I evaluated with the 
best-fit parameters for A^jg and g 2 (A 2 ) (solid blue line). For the sake of comparison, the four-loop perturbative prediction 
evaluated for the same is also plotted with a dashed blue line. 



This is an additional strong constraint, involving the overall bare parameters 



with a^. ert evaluated for the same A^jg 
z(fi 2 , a^ 1 ) , that can be only observed if lattice artefacts are properly treated and Eq. (11) indeed accounts for 
nonperturbative corrections such that z g i (z g h) corresponds 5 to the perturbative bare dressing function for the gluon 
(ghost). Then, the perturbative Taylor coupling can be computed from Eq. (24 |'s r.h.s. at /i = 10 GeV, with the 
results reported in tab. |n] and we thus obtain: a\ 



pert 



0.222(8) for /3 = 1.95 and 

0.211(12) for j3 = 1.90. These values can be compared with a^ crt = 0.214, obtained in perturbation for 



1 ,P crt 



0.216(2) for (3 = 2.10, a 



pert 
T 



A m = 311 MeV 



or with a^ ert 



0.216 for A 



checking of the constraint given by Eq. (24 1. 



MS 



320 MeV. This comparison clearly provides with a pretty positive 



V. CONCLUSIONS 



The gluon and ghost propagators, computed from lattice QCD simulations with two light and two heavy dynamical 
quark flavours, have been successfully described, for momenta above 4.5 GeV, with running formulae including four- 
loop perturbative corrections and a nonperturbative OPE power contribution led by the only dimension-two gluon 
condensate in Landau gauge, g 2 (A 2 ). The OPE formulae including only this leading nonperturbative correction 
fails to describe properly the lattice data for the propagators at momenta below 4.5 GeV, where next-to-leading 
corrections appear to be required. The contribution from the nonperturbative correction to the running is given by 
the Wilson coefficients for the local operator A 2 in the OPE expansions of both gluon and ghost two-point functions. 
After defining the appropriate renormalization scheme for the two-point function and the local operator, the Wilson 
coefficient is also known at the four-loop order. As Aj^rg- for A//=2+l+l is well known from r decays or can be 
consistently obtained from the world average value for the Strong coupling at Zq mass scale, the remaining lattice 
artefacts can be removed and the nonperturbative contribution isolated from data on the same footing. This allows 
for a precise and positive test of the running due to the Wilson coefficients from data, the gluon condensate value 
being the same for both gluon and ghost two-point functions. The main results of the paper appear thus sketched in 
Tab. where the universality for the condensate is checked, and in Fig. [3j where we scrutinize the gluon and ghost 



5 If so, Eq. iL-lli is nothing but the perturbative definition for the Taylor coupling. 
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propagators data to isolate the contribution from the Wilson coefficient to their running with momenta. 

We finally tested that cut directly derived from the bare ghost and gluon dressing functions is consistent with the fits 
on both the renormalized dressing functions; and also compatible A^ obtained from r decays (as shown in ref. |14j). 
This confirms in a non trivial way the validity of our estimates for the lattice artifacts and the nonperturbative 
correction, Eq. (11). 

In summary, we found unequivocal deviations for ghost and gluon lattice propagators with respect to their four-loop 
perturbative prediction that have been consistently and accurately accommodated within the nonperturbative OPE 
approach. 



Appendix A: The matching criterion 



The criterion to determine the best matching for the artefacts-free rescaled dressing functions, simulated at different 
/3 parameters, is based on the minimization of the area between polynomial fits of data. To this purpose, we fit the 
data with Legendre polynomials, 



Gtf^fJ})) = XXflfl 2 / _7 -1 = Ml), (Al) 

^yiiyjJ ) ) . y yraax ymin 



where Pi(x) is the i-th order Legendre polynomial that, being defined within the interval (—1, 1), appears in Eq. ( Al ) 
written in terms of an argument that makes it to range from g min to g max , the lower and upper bounds for the fitting 
window. The coefficient w 4 (/3) corresponds to the weight for the i-th order polynomial. Therefore, we will have: 



n 



For our three lattice data sets, we will take /3q = 2.10 and will fix n + 1 = 10 (the number of terms for the Legendre 
polynomials), g m ; n = 2 GeV and <z max will be the minimum of the largest momenta for each two data sets being 
matched. Then, we finally define 

A p = A(2.10, 1.95) +A(2. 10, 1.90) (A3) 

as the parameter to be minimized. 
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